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1.  INTRODUCTION 


The  conventional  method  of  moments  is  not  suitable  for  electrically  large  bodies  due  to  the 
prohibitively  large  storage  requirements  in  the  computer.  However,  an  alternative  to  directly 
inverting  the  MoM  impedance  matrix  is  to  recast  the  problem  in  a  manner  such  that  iterative 
schemes  can  be  employed  to  solve  the  problem  of  radar  scattering  from  electrically  large  bodies.  In 
this  report,  we  examine  a  number  of  iterative  techniques  that  appear  to  be  useful  for  this  purpose. 

Although  various  types  of  iterative  schemes  for  matirx  equations,  e.g.,  the  conjugate  gradient 
method,  have  been  developed  in  the  past,  they  were  usually  optimized  to  solve  the  problem  for  a 
single  illuminating  field  (right-hand  side)  and,  therefore,  are  not  too  well-suited  for  handling 
multiple  incidences.  However,  since  most  of  the  applications  require  the  computation  of  the  radar 
cross  section  for  a  wide  range  of  incident  angles,  it  is  important  to  be  able  to  treat  the  multiple 
incidence  case  in  an  efficient  manner.  To  this  end,  we  suggest  some  variational-iteration  schemes 
based  on  the  use  of  prechosen,  entire  domain  basis  functions  that  not  only  enable  one  to  treat  larger 
bodies  but  allow  the  handling  of  multiple  incident  angles  as  well.  The  specific  form  of  the  basis 
functions  chosen  as  an  example  in  this  chapter  is  the  traveling-wave  type.  The  advantages  of  using 
this  type  of  basis  functions  are  that  their  functional  form  is  simple,  known  in  advance  and 
independent  of  the  incident  field.  In  this  report,  the  problem  of  TM  scattering  from  edge-loaded 
strips  is  investigated  in  detail  to  illustrate  the  use  of  different  solution  techniques  mentioned  above. 
Extensive  numerical  results  are  presented  and  some  conclusions  are  drawn  on  the  basis  of  these 
results.  The  case  of  TE  scattering  from  edge-loaded  strips  is  also  discussed  briefly. 


2 .  ITERATIVE  METHODS 


vj 


2. 1  Formulation  for  TM  Scattering 

The  geometry  of  the  edge-loaded  strip  and  the  coordinate  system  used  is  shown  in  Figure  2.1 
for  TM  scattering.  On  the  strip,  the  incident  plane  wave  is  given  by 


E‘n(x)  =  z  e^k°xcos9° 


where  Go  is  the  incident  angle  measured  from  the  positive  x-axis.  The  general  three-dimensional 
integral  equation  for  a  thin  resistive  scatterer  has  the  form 


Etm(r)  =  R(r)  J(r)  -  E^r)  r  on  the  surface  S 


where 


Es  =  -jcoA  —  V<b 


A(r)  =  p.f  J(r')G(r,r')  dr’ 

J  s 

0(r)  =  r--— ■  f  V'-J(r')G(r,r')  dr' 
jtoe  js 


where  S  denotes  the  surface  of  the  scatterer.  Specialized  to  an  infinitely  thin  resistive  strip  with  the 
incident  electric  Field  in  TM  polarization,  it  becomes 


Em(x)  =  R(x)  J(x)  -  Es(x)  on  the  strip 


where 


WTFfw  w ~ 


.  *»« 


S 


Es(x) 


=  -jM  J 

•-a 


J(x')  G(x,x')  dx' 


G(x,x')  =  ■^■Ho2)(k0lx-x'l) 


G(a)  = 


W  kJ-a2 


and  ko  denotes  the  free-space  wave  number,  rj  is  the  intrinsic  impedance  of  the  free  space,  and  a 
is  the  Fourier  transform  (spectral)  variable.  In  the  above  equations,  G(x,x')  is  the  two- 
dimensional  Green's  function,  Es(x)  denotes  the  scattered  field,  and  the  quantities  in  italics  indicate 
the  Fourier  transforms  with  a  as  the  transform  variable.  Equation  (2.6)  is  the  electric  field  integral 
equation  (EFIE)  for  the  problem  of  TM  scattering  from  a  resistively  edge-loaded  strip.  The 
convolution  integral  for  the  scattered  field,  given  in  Eq.  (2.7),  can  be  efficiently  computed  using 
the  Fast  Fourier  Transform  (FFT).  In  the  following  sections,  various  ways  of  solving  the  EFIE  in 
Eq.  (2.6)  are  investigated  in  detail. 


2.2  Conjugate  Gradient  Method  (CGM) 

The  conjugate  gradient  method  is  based  on  the  conjugate  direction  method  and  includes  an 
algorithm  for  the  generation  of  mutually  conjugate  p-functions  [1].  It  has  been  applied  to  different 
types  of  electromagnetic  scattering  problems  in  the  past  [2],  [3],  [4],  [5].  In  general,  the  conjugate 
gradient  method  converges  quite  rapidly  for  a  matrix  whose  eigenvalue  spectrum  is  such  that  a 
large  number  of  its  eigenvalues  are  bunched  closely  together.  This  usually  happens  for  diagonally 
dominant  matrices  [6].  For  the  particular  problem  of  TM  scattering  from  an  edge-loaded  strip, 
CGM  works  well  for  any  single  incidence  and  gives  accurate  results  in  a  few  iterations.  However, 
CGM  is  not  the  most  suitable  approach  when  one  desires  to  solve  the  problem  for  many  incident 
angles  simultaneously,  because  one  must  start  anew  and  generate  an  entire  sequence  of  expansion 
functions  (gradient  vectors)  each  time  a  new  incident  field  is  inserted  into  the  integral  equation  such 


,v;v 


VMV 


as  Eq.  (2.6).  Since  CGM  generates  a  sequence  of  the  orthogonal  gradients  which  span  the  entire 
space,  one  might  conjecture  that  a  set  of  gradient  vectors  generated  for  one  incident  angle  could  be 
used  for  all  of  the  other  incident  fields  and,  consequently,  the  generation  of  a  new  set  of  gradient 
vectors  would  not  be  required  each  time  the  angle  of  illumination  is  changed.  In  practice, 
however,  the  numerically-generated  expansion  functions  quickly  lose  their  orthogonality  because 
of  computer  round-off  error  and  it  is  totally  impractical  to  use  the  same  functions  even  for 
moderately  separated  angles  of  incidence.  Since  CGM  is  numerically  rigorous  and  barring  a 
build-up  of  round-off  errors,  it  generates  the  same  solution  as  would  be  obtained  via  the  matrix 
method  in  the  context  of  MoM,  it  is  used  to  check  the  convergence  of  other  iterative  algorithms 
presented  below.  However,  a  detailed  study  of  the  application  of  CGM  and  its  variations  to 
different  types  of  electromagnetic  field  problems  is  beyond  the  scope  of  this  report  and  the  reader  is 
referred  to  [7],  [8]  for  additional  information. 

2.3  Variational  Iteration 

In  this  section,  we  describe  iterative  algorithms  for  simultaneous  treatment  of  multiple  incident 
angles.  We  begin  with  a  set  of  prechosen  expansion  functions  which  give  equal  weight  to  all 
incident  angles.  Next,  we  devise  an  iterative  algorithm  for  minimizing  an  error  functional  which  is 
a  measure  of  the  satisfaction  of  the  boundary  condition  on  the  strip.  Prechoosing  the  direction 
vectors  (expansion  functions)  obviates  the  need  for  generating  them  repeatedly  for  different 
incident  angles,  as  is  done  in  CGM.  In  this  section,  the  expansion  functions  are  chosen  to  be  the 
plane-wave  type  although  one  could  just  as  easily  choose  them  differently,  if  one  so  desired.  After 
the  iteration  is  carried  out,  the  current  is  expressed  in  the  form 

N 

J(x)  =  Jc(x)  +  £  Tin  e)koxcos9n  (2. 1 0) 

n=l 


,  A  *\  «\  f  ^ 


,■»  .V.V. 
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where  ko  is  the  free-space  wave  number,  0n  is  the  angle  measured  from  the  positive  x-axis,  J0(x) 
is  an  initial  guess,  and  r)n  is  the  complex  weighting  coefficient  generated  in  a  manner  to  be 
explained  below.  Once  these  weighting  coefficients  are  determined,  we  have  the  traveling-wave 
representation  of  the  induced  current  on  the  strip.  There  are  two  advantages  to  using  plane-wave- 
tvpe  basis  functions.  First,  their  functional  form  is  very  simple  and,  consequently,  they  can  be 
computed  easily  and  inexpensively.  Second,  since  the  current  is  already  expressed  in  terms  of 
traveling  waves,  it  is  a  straightforward  matter  to  identify  which  traveling-wave  component  has  the 
dominant  contribution  to  the  far  field  for  a  given  angle  of  observation.  To  illustrate  this  point, 
suppose  we  start  with  zero  as  the  initial  guess  of  the  current  on  the  strip.  Then  we  have 

N 

J(x)  =  Spn^  (2.11J 

n=l 


The  scattered  field  in  the  far  zone  can  be  found  by  using  the  asymptotic  expression  for  the  Hankel 
function  [9]  and  can  be  written  as 


In 


sin[koa(cos9+cos0n)] 

koa(cos0+cos9n) 


(2.12) 


Equation  (2  12)  shows  that  the  far  field  is  a  combination  of  sine  functions,  i.e.,  sinx/x  type 
functions.  These  functions  peak  at  different  angles,  given  by  (n-0o),  where  the  argument  of  the 
sine  function  becomes  zero  Note  that,  in  a  method  such  as  CGM,  it  may  be  quite  difficult  to 
extract  different  traveling- wave  components  from  the  final  induced  current  obtained.  We  also 
observe  that  when  the  basis  functions  are  prechosen,  i.e.,  they  are  not  generated  via  the  iteration 


algorithm,  the  same  set  of  basis  functions  can  be  employed  for  an  arbitrary  incident  angle. 
Although  this  kind  of  generality  in  the  method  is  a  very  desirable  characteristic,  it  should  be 
pointed  out  that  the  accuracy  we  obtain  for  each  incident  angle  may  not  be  as  good  as  that  obtained 
by  CGM  or  any  other  single-incidence  oriented  algorithm.  However,  as  has  already  been  noted 
before,  these  single-incidence  oriented  algorithms  generate  gradient  vectors  that  are  optimal  only 
for  a  particular  incident  angle,  and,  consequently,  are  ill-suited  for  the  treatment  of  multiple 
incident  angles.  In  contrast,  the  use  of  an  iterative  algorithm  with  prechosen  basis  functions 
allows  one  to  treat  multiple  incident  angles  simultaneously  with  the  same  set  of  basis  functions. 
The  conventional  method  of  moments  may  also  be  employed  effectively  if  the  matrix  size  is  not  too 
large.  The  purpose  of  developing  an  iteiative  procedure  is  to  circumvent  the  problem  of  large 
matrices  which  cannot  be  conveniently  treated  by  the  direct  inversion  method  and  to  exploit  any 
symmetry  or  sparseness  that  arises  due  to  the  nature  of  the  problem  under  consideration. 

2.3. 1  Minimization  of  the  boundary  condition  error 

In  this  section,  an  algorithm  for  variational  iteration  with  the  second  minimization  step  (SMS), 
minimizing  I  y  -  Lxa  |2,  is  derived,  where  y  is  the  incident  field,  xa  is  an  approximate  solution,  and 
L  is  the  integral  operator.  The  derivation  of  SMS  for  this  particular  type  of  functional  has  been 
given  in  [5]  and  will  not  be  repeated  here.  The  resultant  expressions  are  summarized  below. 

Suppose  we  are  given  a  set  of  basis  functions  gj  through  gN  and  an  integral  equation  of  the 
form 

L  xex  =  y  (2.13) 

where  xex  is  the  unknown  to  be  determined  As  the  result  of  minimizing  the  functional  on  the 
boundary  condition  error,  the  weighting  coefficients  for  these  basis  functions  are  given  by 


v  *•  n-l»  *n  " 

<  fn.  fn  > 


(2.14) 


where 


Fn-1  =  y  -  L 


fn  =  L  §n 


§n  §n  §n-l 


(2.15) 


fn  =  fn  -  fn-l 


<fn-l.fn> 


<fn-l.fn-l> 


where  gn's  are  the  predetermined  basis  functions,  while  gn’s  with  the  bar  on  top  denote  the 
modified  gn's  obtained  via  the  second  minimization  step  (SMS).  The  inner  product  is  defined  as 


f.  g  >  =  f  f'g* 

Jo 


(2.16) 


where  D  is  the  domain  of  the  given  problem  and  denotes  the  complex  conjugate. 


2.3.2  Minimization  of  the  err^r  on  the  induced  cunrent 


The  convergence  of  the  iterative  algorithm  minimizing  the  boundary  condition  error,  as 
described  in  Section  2.3.1,  does  not  always  guarantee  the  convergence  of  the  approximate  solution 


5] 


to  the  true  one.  It  has  been  noted  in  the  past  that  situations  exist  where,  although  the  boundary 
condition  error  is  reduced  below  the  prescribed  level,  the  computed  solution  is  not  close  to  the  true 
solution.  In  this  section  an  algorithm  for  variational  iteration  with  the  second  minimization  step, 
which  directly  minimizes  the  error  on  the  solution  itself,  i.e.,  I  xex  -  *a  *s  derived,  where  xex 
denotes  the  exact  solution  to  the  problem  under  consideration  and  xa  an  approximate  solution 
which  is  generated  in  a  suitable  manner.  The  minimization  of  the  error  on  the  current  is  achieved 
by  utilizing  the  adjoint  operator,  L^.  Since  we  do  not  know  the  exact  solution,  xex  .  to  the 
problem,  we  cannot  explicitly  calculate  the  quantity  1  xex  -  xa  —  However,  it  is  still  possible  to 
minimize  the  functional  I  xex  -  xa  |2  .  Since  the  derivation  for  SMS  minimizing  this  functional  is 
quite  different  from  that  given  in  Section  2.3.1,  where  the  boundary  condition  error  is  minimized, 
we  show  below  the  steps  involved  in  detail. 

Suppose  that  we  have  already  iterated  up  to  the  (n-l)th  step.  Then  we  have 


^n-1  Xex  ^n-l 


(2.17) 


x  ^n-l >  ^“n-l  '>  x  Xex  Xex  ^n-l 


(2.18) 


where  Fn  denotes  the  residual  on  the  current  and  ERRn  the  value  of  the  error  evaluated  at  the  nfb 
step.  For  the  problem  of  TM  scattering  from  the  edge-loaded  strip,  the  adjoint  operator  is  known 
and  is  given  by 


La(  )  =  R(x)  ( 


•)-f  (  )  G*(x, 
Jd 


x1)  dx' 


(2.19) 


where  means  the  complex  conjugate  of  the  quantity.  Let  the  expansion  functions  'f'n  be  the 
plane-wave  functions,  i.e.. 


Wv< 


rij" 
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Notice  that  <  xn- 1,  gn  >  is  known  since  xn- 1  and  gn  are  known.  Also 


<  Xex»  gn  >  =  <  Xex-  L"Vn  >  =  <  Lxex,  \|/n  >  =  <  y,  Vn  > 


Therefore,  we  can  write 


<  F„-I.  gn  >  _  <  y.  Vn  >  -  <  *n-l’  L-Vn  > 

<  gn’  gn  >  <  LVn,  L\n  > 


(2.26) 


(2.27) 


All  the  terms  appearing  in  the  above  expression  for  T}n  are  known.  Having  found  the  weighting 
coefficient  rjn,  we  can  rewrite  the  error  given  in  Eq.  (2.23)  in  the  following  manner.  For  the 
convenience  of  notation,  let 

An  =  <  Fn-l.  gn  >  =  <  y  -  Ut-l-  Vn  >  (2.28) 

®n  —  ^  gn’  gn  ^ 

where  the  current  expansion  functions,  gn's,  are  given  by  Eq.  (2.21).  Then  the  error  on  the 
current  is  given  by 


ERRn  =  <Fn_1,Fn_1>- 


An_ 

B„ 


(2.29) 


The  second  minimization  step  aims  at  maximizing  the  second  term  in  Eq.  (2.29),  thereby  obtaining 
the  greatest  reduction  possible  in  the  error.  One  way  of  achieving  such  maximization  is  to  define  a 
new  gradient 
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Sn  Bn  ^nSn-1 


(2.30) 


where  is  the  unknown  weighting  coefficient  to  be  determined  as  the  result  of  the  second 
minimization  step.  Since  we  can  write  the  residual  as 


—  _  _  c  _  **  Fn-1»  8n  > 

Fn  —  *ex  —  Fn— 1  ^InSn  —  Fn-1  §n 

<  En  >  gn  > 


we  have 


(2.31) 


r-  _  r,  _  Fn-1»  §n  ^ 

<  Fn>  §n  '>  ~  ^  Fn-1>  §n  ^  —  <"  Bn’  §n  ^ 

<  8n»  8n  ** 


(2.32) 


=  <  Fn-1>  Sn  >  “  <  Fn-l.  Sn  >  =  0 


Notice  that  the  numerator  of  the  second  term  in  Eq.  (2.29)  remains  the  same,  i.e.. 


^  Fn-1*  8n  ^nSn-1  ^ 


^  Fn-1»  Bn  ■>  Fn-1>  En-1  ^  <  Fn-1»  §n  ^  ~ 


(2.33) 


The  modified  denominator  is  given  by 


®n  ^  8n’  Sn  ■> 


(2.34) 


^  §n>  §n  ^  ^  §n— 1>  §n  ^  ^  Bn’  §n-l  ^  ^  8n-l*  §n-l  ^ 


Since  the  numerator  remains  the  same,  the  maximization  of  the  second  term  in  Eq.  (2.29)  is 
achieved  by  minimizing  the  denominator.  Upon  carrying  out  this  minimization,  we  obtain  the 
following  expression  for  the  unknown  weighting  coefficient  : 


^  Sn  Sn-1  *> 
<'  §n—  1’  Sn-1  ^ 


Inserting  Eq.  (2.35)  into  Eq.  (2.34),  we  finally  obtain 


Bn=Bn- 


<  Sn»  Sn-l^ 

B„-t 


(2.35) 


(2.36) 


This  completes  the  implementation  of  the  variational-iteration  procedure  for  minimizing  the  residual 
on  the  current.  As  has  been  noted  before,  the  error  on  the  current  itself  cannot  be  computed 
explicitly  since  we  do  not  have  access  to  the  true  solution.  Thus,  to  check  the  convergence  of  the 
iteration  procedure,  the  functional  form  of  the  approximate  current  at  each  step  of  the  iteration  is 
substituted  into  the  expression  for  the  error  on  the  boundary  condition,  and  this  error  is  used  as  a 
measure  of  the  accuracy  of  the  solution  for  the  current. 


2.3.3  Numerical  results 

Figures  2.2,  2.3  and  2.4  show  the  bistatic  radar  cross  sections  obtained  for  three  different 
incident  angles  by  the  variational-iteration  methods  described  above  and  they  are  compared  to  those 
computed  from  CGM.  The  strip  width  is  2  X  and  is  edge-loaded  with  the  maximum  resistivity  of 
200  ohms  at  the  edges  of  the  strip  and  tapered  in  a  linear  fashion  to  zero  resistivity  at  0.41  X  from 
each  edge.  In  the  variational  iteration,  the  same  set  of  nine  plane-wave-type  basis  functions  is  used 
as  the  prechosen  gradient  vectors  for  all  three  incident  angles.  As  can  be  seen  in  the  figures,  two 
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variational  iteration  schemes  described  in  Sections  2.3.1  and  2.3.2  seem  to  yield  about  the  same 
results  for  the  particular  problem  of  TM  scattering  from  edge-loaded  strips,  in  terms  of  accuracy 
and  convergence. 

2.4  Preconditioning  in  the  Spectral  Domain 
Preconditioning  of  a  matrix  equation 

A  x  =  b  (2.37) 

where  A  is  the  matrix  operator,  b  the  excitation  vector,  and  x  the  unknown,  usually  consists  of 
finding  a  pair  of  matrices  P  and  Q  and  transforming  the  original  equation  into  an  equivalent  pair  of 
equations 


PAQy  =  Pb 

(2.38) 

X 

II 

€ 

(2.39) 

If  the  matrices  P  and  Q  are  chosen  such  that  the  condition  number  of  PAQ  is  smaller  than  that  of  A, 
an  iterative  algorithm  such  as  CGM  will  generally  converge  faster  on  the  preconditioned  system 
than  on  the  original  system  [10],  [11],  In  this  section,  a  simple  way  of  preconditioning  an  integral 
equation  in  the  spectral  domain  is  introduced.  The  advantage  of  doing  the  preconditioning  in  the 
spectral  domain  is  that  the  operation  involved  is  a  simple  algebraic  multiplication  or  division  rather 
than  matrix  multiplication  or  matrix  inversion.  The  integral  equation  in  the  space  domain  is  given 
by 


R(x) 


K,)-f 


J(x') 


Gp(x,x’)  dx’  =  y(x) 


(2.40) 


m 


on  the  strip,  where  Gp  is  the  modified  Green's  function,  which  is  simply  the  original  Green’s 
function  multiplied  by  an  appropriate  scale  factor.  Let  H  be  a  preconditioning  operator  in  the 
spectral  domain,  Fourier  transforming  the  above  integral  equation,  multiplying  both  sides  of  the 
transformed  equation  by  H  and  inverse  Fourier  transforming  the  resultant  equation  yield  the 
preconditioned  version  of  Eq.  (2.40).  The  preconditioned  operator  equation  can  be  written  as 

0F-1  <  H  F[0R(x)  J(x)]  >  -  0F_1  <  H  F[0F_1C/Gp)]  >  =  0F"1  <  //F[0y]  >  (2.41) 

where  0  is  a  truncation  factor  which  is  set  to  1  if  on  the  strip  and  0  otherwise.  F  denotes  the 
forward  Founer  transform  and  F*  1  the  inverse  Fourier  transform.  The  idea  is  to  choose  H  such 
that  the  preconditioned  operator  will  approximate  an  identity  operator.  If  the  strip  is  perfectly 
conducting,  then  choosing  the  preconditioning  operator  H  to  be  the  inverse  of  Gp  guarantees  the 
preconditioned  operator  to  be  an  approximate  identity  operator,  and  an  iterative  algorithm  applied  to 
this  preconditioned  equation  converges  very  rapidly.  However,  if  the  strip  is  edge  loaded,  it  can 
be  shown  that  there  is  no  general  way  of  finding  the  operator  H  that  will  lead  to  an  approximate 
identity  operator  after  preconditioning.  Also,  it  is  noted  that  when  the  operator  is  preconditioned  in 
the  manner  described  in  Eq.  (2.41),  the  use  of  prechosen  gradient  vectors,  such  as  plane-wave- 
type  functions,  yields  poor  results.  To  have  reasonable  results  with  this  type  of  preconditioning,  it 
is  necessary  to  choose  the  gradients  that  are  closely  related  to  the  residual  Fn-l,  which,  in  turn,  is 
intimately  related  to  the  incident  angle.  Therefore,  as  in  the  case  of  CGM,  this  type  of 
preconditioning  in  the  spectral  domain  may  not  be  suitable  to  treat  multiple  incident  angles. 

Some  of  the  results  for  a  4.5  X  wide  PEC  strip  obtained  by  the  preconditioning  technique 
described  above  are  shown  in  Figures  2.5,  2.6,  2.7  and  2.8.  Preconditioning  is  done  with  the 
inverse  of  Gp  as  the  preconditioning  operator.  The  current  expansion  functions  are  chosen  to  be 
the  residuals,  Fn- 1  ’s.  When  the  strip  is  perfecdy  conducting,  preconditioning  in  this  manner  gives 
faster  convergence  than  CGM.  However,  as  Rmax  increases,  the  number  of  iterations  needed  to 
achieve  the  same  order  of  accuracy  as  that  given  by  CGM  also  increases. 
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Magnitude  of  the  induced  current  for  normal  incidence. 
Solid  line  —  CGM. 

Dotted  line  —  Preconditioning. 


Bistatic  RCS  for  normal  incidence 
Solid  line  —  CGM. 

Dotted  line  —  Preconditioning. 
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2.5  Spectral  Iterative  Technique  (SIT) 


2.5.1  Conventional  SIT 


The  spectral  iterative  technique  has  been  applied  to  electromagnetic  problems  by  a  number  of 


authors  [12],  [13],  [14],  [15],  [16],  The  conventional  SIT  algorithm  can  be  derived  for  the  edge- 


loaded  strip  case  as  follows.  The  integral  equation  in  Eq.  (2.40)  can  be  rewritten  as 


0R(x)  J(x)  -  0  f  J(x')  Gp(x,x')  dx1  =  0y(x)  for  all  x 

J  — oo 


(2.42) 


where  0  is  the  truncation  factor  mentioned  before.  Rearranging  the  above  equation,  we  get 


A  -OO  —oo 

0R(x)  J(x)  -  0y(x)  +  0  J  J(x')  Gp(x,x’)  dx'  =  J  J(x')Gp(x,x')  dx' 


(2.43) 


where  9  is  the  complementary  truncation  factor  which  is  0  on  the  strip  and  1  otherwise.  Since 


f  J(x')  Gp(x,  x’)  dx’  =  r  \JGV  ) 

* — oo 


(2.44) 


we  can  write  Eq.  (2.43),  after  Fourier  transforming  it,  as 


JGp  =  F<0R(x)  J(x)  -  0y(x)  +  0F_1C/Gp)> 


(2.45) 


From  this  equation  we  can  get  the  iteration  scheme 


4+1  =  c;1  F<  0R(x)  Jn(x)  -  0y(x)  +  0  r1(4GJ  > 


(2.46) 


vyv.y.v.v,.  * . 


nWLVT.V'. 
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This  is  the  conventional  SIT  equation  solving  for  the  current  induced  on  the  strip.  The  success  of 
this  algorithm  depends  upon  the  assumption  that  Jn+1»  computed  from  Eq.  (2.46),  represents  a 
better  approximation  to  the  true  solution  than  the  previous  iterate  Jn.  Notice  that,  if  the  strip  is  not 
perfectly  conducting,  then  the  inverse  of  G p  does  not  approximate  too  well  the  true  inverse 
operator  of  the  integral  equation  given  in  Eq.  (2.42)  and  the  solution  does  not  necessarily  improve 
with  iteration.  In  fact,  when  R(x)  exceeds  a  certain  maximum  value,  the  algorithm  actually 
becomes  unstable  or  divergent.  However,  for  small  values  of  R(x),  SIT  appears  to  work  quite 
well  because,  as  the  resistivity  on  the  strip  is  reduced  to  zero  and  the  strip  becomes  perfectly 
conducting,  the  inverse  of  Gp  does  indeed  approximate  quite  closely  the  true  inverse  of  the  original 
operator. 


2.5.2  Modified  SIT  :  SIT  combined  with  SMS 

The  spectral  iteration  with  the  second  minimization  step  is  an  attempt  to  improve  the  original 
SIT  algorithm  such  that  the  procedure  does  not  diverge  even  at  its  worst.  Of  course,  this  does  not 
mean  the  algorithm  will  converge  quickly.  However,  although  the  convergence  may  be  slow  in 
some  cases,  this  combination  guarantees  that  the  algorithm  will  not  diverge  as  in  the  case  of  the 
original  SIT.  The  basic  idea  is  to  use  the  SIT  equation  to  generate  the  gn's,  the  current  expansion 
functions  instead  of  generating  the  Jn's,  the  updated  currents  themselves.  The  equation  for 
cal  mg  the  gn's  is  given  by 


new  gn  =  G  1  F<  9R(x)  gn(x)  -  9Fn_;(x)  +  9F  !(gn  Gp)  > 


(2.47) 


where 


Fn  =  Y  “  L  xn 


(2.48) 


SSKSSig 
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With  the  gn's  obtained  from  Eq.  (2.47),  we  proceed  to  SMS  to  compute  the  appropriate  rjn's,  the 
weight  coefficients,  for  these  current  expansion  functions  by  minimizing  the  boundary  condition 
error.  Although  improved  in  the  sense  that  it  does  not  diverge  for  all  the  values  of  resistivity, 
results  show  that  the  convergence  of  this  modified  SIT  becomes  slow  when  R(x)  exceeds  a  certain 
maximum  value. 


2.6  Block  iteration  :  A  Hybrid  Approach  Based  on  a  Combination  of  Iteration  with  MoM 

To  solve  radar  scattering  problems  for  electrically  large  bodies  by  conventional  methods 
requires  large  storage  and  computation  times;  these  are  major  limiting  factors  in  many  cases.  The 
block  iteration  approach  is  a  scheme  which  attempts  to  circumvent  these  difficulties.  The  idea  of 
block  iteration  is  not  new.  In  the  past,  however,  the  conventional  block  iteration  had  convergence 
problems,  i.e..  the  procedure  diverged  in  many  cases.  In  this  section,  it  is  shown  that  the  block 
iteration  can  always  be  made  nondivergent. 

The  block  iteration  is  a  scheme  for  solving  a  system  of  linear  equations  iteratively  by  scanning 
the  matrix  with  a  small  window.  It  requires  the  matrix  to  be  divided  into  small-order  blocks,  or 
windows,  along  the  main  diagonal.  These  windows  are  then  solved  in  sequential  manner,  and  the 
appropriate  coefficients  obtained  by  minimizing  the  boundary  condition  error  are  updated  each  time 
a  window  is  solved.  The  advantage  lies  in  the  fact  that  the  entire  matrix  need  not  be  stored.  At 
each  step  of  the  iteration  procedure,  only  the  window  must  be  stored.  The  scheme  is  explained  by 
the  following  simple  example. 

Suppose  we  have  four  prechosen  basis  functions  gi,  g2,  g3  and  g4.  Given  a  problem  of  the 


L  x  =  y 


(2.49) 


where  L  is  a  linear  operator,  y  is  a  given  excitation  function,  and  x  the  unknown  to  be  determined, 
we  start  the  iteration  procedure  with  an  initial  guess  xD  obtained  in  a  suitable  manner.  Typically 
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this  initial  guess  is  made  by  using,  e.g.,  the  physical  optics  approximation,  or  by  simply  choosing 
it  equal  to  zero.  Next,  we  redefine  the  unknown  as  the  difference  between  the  original  unknown, 
x,  and  the  initial  guess,  Xq,  i.e.. 


Lxi  =  y-Lx0  =  y1 


(2.50) 


where  xi  is  now  the  unknown  to  be  determined,  and  solve  Eq.  (2.50)  with  the  basis  functions  gi 
and  g2  by  using  the  method  of  moments  to  obtain 


Xl^l^gl+Tl^ga 


(2.51) 


where  the  subscripts  denote  the  corresponding  basis  functions  and  the  superscript  the  stage  of 
iteration.  Using  this  xj,  we  go  through  the  usual  SMS  procedure  outlined  in  Section  2.3.1  and 
calculate  a  weight  coefficient  CU)  which  minimizes  the  error  functional  I  yi  -  Lxj  |2.  Using  this 
weight  coefficient,  we  write  an  improved  solution  to  Eq.  (2.50)  as 


x^C^Vi^gi+Ti^) 


(2.52) 


To  proceed  to  the  next  stage  of  iteration,  we  again  redefine  our  unknown,  X2,  to  be  the  difference 
between  the  original  unknown,  x,  and  the  sum  of  Xq  and  xj,  and  write 


L  x2  ~  y  Lx0  L  Xj  =  y  Lx0  —  ^Cqj^Lgi  +  Tj^Lgj)  =  ] 


(2.53) 


By  solving  this  equation  with  g2  and  g3  and  repeating  the  above  procedure,  we  get 


X  _  /-'(2)  /  (2)  (2)  1 

X2  -  L  (T)2  g2  +  T)3  g3) 


(2.54) 
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This  completes  the  second  stage  of  the  iteration.  Finally,  repeating  the  above  procedure  once  more 
and  adding  an  initial  guess  and  those  three  "partial"  solutions  at  each  of  the  three  iteration  stages, 
we  obtain  the  solution  to  the  original  equation 

x  =  x„  +  C™  v/'g,  +  (  c%'”  +  C%f )  g2 

+  (C®T,f  +  C(3>T)?>)g,  +  Cl3,ll?,g4  (2.55) 

which  is  the  desired  solution  obtained  by  the  block  iteration  with  the  basis  functions  gi  through  g4 
and  with  the  window  size  of  two.  In  short,  the  block  iteration  is  a  scheme  for  generating 
appropriate  expansion  functions  for  the  solution,  the  xn's  in  the  above  example,  by  utilizing  the 
method  of  moments  and  then  improving  the  accuracy  of  the  solution  through  the  SMS  procedure 
in  a  sequential  manner.  Each  of  these  xn's  is  a  certain  combination  of  the  basis  functions,  gn's, 
determined  by  the  method  of  moments.  Notice  that  if  the  gn's  are  subdomain  basis  functions,  then 
"uncoupling"  between  the  matrix  elements  may  occur  and  the  block  iteration  achieves  fast 
convergence.  Howvever,  if  the  gn's  are  entire  domain  basis  functions  that  are  not  well  localized, 
as  in  the  case  of  plane-wave-type  basis  functions,  there  is  no  uncoupling  between  the  matrix 
elements,  and  the  iteration  procedure  converges  slower.  Some  results  obtained  by  using  the  block 
iteration  with  the  plane-wave-type  basis  functions  show  that  this  type  of  hybrid  approach,  where 
an  iteration  procedure  is  combined  with  the  method  of  moments,  is  quite  promising.  As  can  be 
expected,  the  solution  improves  as  the  size  of  the  window  is  increased. 

2.7  Initial  Guess 

In  using  an  iterative  scheme,  such  as  block  iteration  explained  in  the  previous  section,  it  is 
important  to  start  the  algorithm  with  a  good  initial  guess  if  a  rapid  convergence  is  to  be  achieved. 
There  may  be  several  different  ways  of  providing  a  reasonably  good  initial  guesses  for  a  particular 
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problem.  In  the  following,  we  show  some  simple  and  inexpensive  ways  to  generate  reasonably 
good  initial  guesses  by  using  the  problem  of  scattering  from  edge-loaded  strips  as  an  example. 

2.7,1  Angle-stepping 

In  the  angle-stepping  scheme,  we  first  generate  a  good  solution  for  one  incident  angle  by  using 
a  numerically  rigorous  method,  such  as  the  method  of  moments  or  CGM,  and  then  use  that 
solution  as  an  initial  guess  for  another  incident  angle  of  interest.  For  example,  a  good  solution  for 
normal  incidence  may  be  generated  and  used  as  an  initial  guess  for,  say,  10  degrees  off  the  normal 
incidence.  The  size  of  the  angular  increments  can  be  varied  to  achieve  the  desired  accuracy  in  the 
solution. 


2.7.2  Physical  optics 

For  resistively-loaded  strips,  the  physical  optics  approximation  implies  the  following  : 


if  R(x)  >  y 


2nxHm(x)  ifR(x)<y 


(2.56) 


where  tj  is  the  intrinsic  impedance  of  the  free  space  and  n  the  unit  normal  vector  to  the  surface  of 
the  strip.  The  physical  optics  approximation  provides  a  fairly  good  initial  guess  if  the  incident  field 
on  the  strip  is  close  to  normal,  and  if  the  strip  is  perfectly  conducting.  However,  when  one  of 
these  conditions  is  not  satisfied,  the  physical  optics  approximation  is  usually  quite  poor. 


2.7.3  Use  of  adjoint  operator 

Another  way  of  generating  an  initial  guess  of  the  edge-loaded  strips  is  to  make  use  of  the 
adjoint  operator.  Given  an  operator  equation 

L  x  =  y  (2.57) 

we  compute  the  initial  guess  to  be 

xc  =  tj0  LA  y  (2.58) 

where  LA  denotes  the  adjoint  operator,  y  the  given  excitation  function,  and  T]0  a  complex  constant 
obtained  as  a  result  of  minimizing  the  boundary  condition  error  at  the  initial  step  of  iteration.  The 
effectiveness  of  this  type  of  initial  guess  depends  on  how  well  the  adjoint  operator  approximates 
the  inverse  operator  for  a  given  angle  of  incidence.  For  scattering  from  the  edge-loaded  strips,  it 
has  been  found  that  this  choice  of  initial  guess  is  somewhat  complementary  to  the  physical  optics 
approximation.  That  is,  this  method  yields  a  relatively  good  initial  guess  if  the  angle  of  incidence 
is  close  to  grazing  and  if  the  strip  is  heavily  edge- loaded,  i.e.,  Rmax  is  very  large.  Figures  2.9  and 
2.10  show  that  the  adjoint  operator  LA  approximates  the  inverse  operator  L*1  quite  closely  for  the 
grazing  incidence  case.  The  dotted  line  represents  the  solution  obtained  by  LAy.  Notice  that  it 
differs  from  the  rigorous  solution  generated  by  CGM  by  only  a  complex  constant.  When  this 
complex  constant  is  generated  by  minimizing  the  boundary  condition  error  I  y  -  Lxa  K  and 
multiplying  by  LAy,  the  resultant  curve  closely  approximates  the  solution  by  CGM  as  shown  in  the 
figures.  The  strip  is  4.5  X  wide,  with  edge  load  width  1  X  from  each  edge,  and  with  Rmax  equal 
to  200  ohms. 

For  all  of  the  three  different  methods  of  generating  an  initial  guess  given  above,  the 
multiplication  of  a  complex  coefficient  obtained  as  a  result  of  minimizing  the  boundary  condition 
error  at  the  initial  step  of  iteration  will,  in  general,  improve  the  initial  guess.  Also,  when  the  restart 
of  the  iteration  is  used,  multiplying  the  solution  by  the  coefficient  obtained  from  minimizing  the 
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boundary  condition  error  at  the  last  stage  of  the  previous  set  of  iterations  will  improve  the 
convergence  of  the  algorithm. 

2.7.4  Numerical  results 

The  following  figures  show  the  effect  of  combining  the  block  iteration  with  different  types  of 
initial  guesses.  The  strip  is  4.5  X  wide,  the  width  of  the  edge  load  is  1  X  from  each  edge,  and 
Rmax  is  equal  to  200  ohms.  Figures  2.1 1  and  2.12  are  obtained  by  block  iteration  using  the  PO 
approximation  as  an  initial  guess  for  the  grazing  incidence  case.  Figures  2.13  and  2.14  are  the 
results  for  the  grazing  incidence  case  obtained  by  stepping  from  the  solution  for  normal  incidence 
generated  by  CGM  to  grazing  incidence  with  an  angular  increment  of  10  degrees.  Notice  the 
siginicant  improvement  achieved  by  stepping  as  compared  to  the  case  where  the  PO  approximation 
is  used  as  an  initial  guess. 

The  size  of  the  angular  increment  and  the  range  of  stepping  are  two  important  factors  in  the 
scheme  of  angle-stepping  as  demonstrated  in  the  following  examples.  The  strip  is  10.5  X  wide  and 
perfectly  conducting.  Figure  2.15  shows  the  solution  for  the  grazing  incidence  case  obtained  by 
MoM  using  twenty  plane-wave-type  basis  functions.  Notice  that  the  MoM  solution  does  not 
approximate  too  well  the  singular  behavior  of  the  current  at  the  edges.  Figure  2.16  shows  the 
solution  generated  by  block  iteration  by  stepping  from  the  CGM  solution  for  the  normal  incidence 
case  to  grazing  incidence  with  angular  increments  of  10  degrees.  The  number  of  the  basis 
functions  is  kept  constant  at  twenty.  Notice  the  improvement  achieved  both  in  the  current  and  in 
the  far  field.  Figure  2.17  shows  the  results  obtained  by  beginning  with  the  CGM  solution  for  20 
degrees  off  grazing  and  using  angular  increments  of  5  degrees  to  step  to  the  grazing  incidence  case. 
By  reducing  the  angular  increment  and  the  range  of  stepping,  a  better  improvement  is  achieved 
both  in  the  current  and  in  the  far  field.  Finally,  as  shown  in  Figure  2.18,  the  agreement  is  very 
good  when  we  begin  with  the  170  degree  CGM  solution  and  step  to  the  grazing  incidence  with 
angular  increments  of  2  degrees. 
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Magnitude  of  the  induced  current  for  grazing  incidence. 

Solid  line  —  CGM. 

Dotted  line  —  Block  iteration  using  1 5  plane-wave-type  functions 
with  window  size  5  and  stepping  from  normal  to 
grazing  incidence  with  angular  increment  of  10  degrees. 


Figure  2.14  Bistatic  RCS  for  grazing  incidence. 

Solid  line  —  CGM. 

Dotted  line  —  Block  iteration  using  15  plane-wave-type  functions 
with  window  size  5  and  stepping  from  normal  to 
grazing  incidence  with  angular  increment  of  10  degrees. 
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Figure  2.16  Top:  Magnitude  of  the  induced  current  for  grazing  incidence. 

Bottom:  Bistatic  RCS  for  grazing  incidence. 

Solid  line  —  CGM. 

Dotted  line  --  Block  iteration  using  20  plane- wave-type  functions 

stepping  from  normal  to  grazing  incidence  with  angular 
increment  of  10  degrees. 
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Figure  2.18  Top:  Magnitude  of  the  induced  current  for  grazing  incidence. 
Bottom:  Bistatic  RCS  for  grazing  incidence. 

Solid  line  —  CGM. 

Dotted  line  --  Block  iteration  using  20  plane-wave-type  functions 
stepping  from  170  degrees  to  grazing  incidence  with 
angular  increment  of  2  degrees. 


The  situation  is  similar  for  the  case  of  edge-loaded  strips,  i.e.,  we  can  obtain  better  results  by 
stepping  in  angle,  just  as  in  the  case  of  PEC  strips.  We  present  the  results  for  a  10.5  X  wide  strip 
with  the  edge  load  width  1  X  from  each  edge  and  Rmax  equal  to  200  ohms.  Figures  2. 19  and  2.20 
show  the  MoM  solution  and  the  solution  derived  by  stepping,  respectively,  with  the  latter  obtained 
by  starting  with  the  normal  incidence  CGM  solution  and  incrementally  stepping  all  the  way  to  the 
grazing  incidence.  It  is  especially  worth  noting  that  the  disagreement  between  the  current  solutions 
at  the  edges  of  the  strip,  generated  by  CGM  and  MoM,  is  almost  entirely  eliminated  when  stepping 
is  used. 

2.8  Formulation  for  TE  Scattering 

The  geometry  of  an  edge-loaded  strip  and  the  coordinate  system  used  for  the  TE  case  is  shown 
in  Figure  2.21.  The  incident  fields  are  given  by 

Jko(xcos0o  +  ysin0o) 

Hm(x,y)  =  z  - - - -  (2.59) 

E^x.y)  =  (  x  sin0o  -  y  cos90  )  e^*0*.*  ^ 

where  0O  is  the  incident  angle  measured  from  the  positive  x-axis  and  r\  the  free-space  wave 
impedance.  For  the  TE  polarization,  only  the  x-component  of  the  current  on  the  strip  is  nonzero. 
When  the  general  three-dimensional  integral  equation  for  a  resistive  scatterer  is  specialized  to  an 
infinitely  thin  resistive  strip  case,  it  becomes 

Ein(x)  =  R(x)  Jx(x)  -  Ej(x)  (2.60) 


on  the  strip,  where  the  incident  field  is  given  by 
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Figure  2.19  Top:  Magnitude  of  the  induced  current  for  grazing  incidence. 
Bottom:  Bistatic  RCS  for  grazing  incidence. 

Solid  line  —  CGM. 

Dotted  line  --  MoM  using  20  plane-wave-type  functions. 
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Figure  2.20  Top:  Magnitude  of  the  induced  current  for  grazing  incidence. 
Bottom:  Bistatic  RCS  for  grazing  incidence. 

Solid  line  —  CGM. 

Dotted  line  -  Block  iteration  using  20  plane-wave-type  functions 
stepping  from  normal  to  grazing  with  angular 
increment  of  10  degrees. 
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Ein(x)  =  sin90  ejk«xcos0° 


(2.61) 


The  scattered  field  can  be  written  as 


E^(x)  =  ~  +  kj)f  Jx(x')  G(x,x')  dx’ 

V  dx  J 


(2.62) 


where  is  the  free-space  wave  number  and  G(x,x')  is  the  two-dimensional  Green's  function 


given  by 


G(x,x')  =  H^(k0lx-x'l) 


(2.63) 


Fourier  transforming  both  sides  of  Eq.  (2.63),  we  have 


£^(a)  =  jtue  (k°  ”  °-2)  G(a)  = 


rkJ-a2 


JJa) 


(2.64) 


where  a  is  a  spectral  variable.  By  inverse  Fourier  transforming  the  right-hand  side  of  Eq.  (2.64), 
we  obtain  the  scattered  field  in  the  space  domain.  All  the  techniques  that  are  used  to  solve  for  TM 
scattering  also  apply  to  the  TE  case.  In  this  section,  only  the  results  obtained  from  the  conjugate 
gradient  method  and  the  block  iteration  using  traveling  waves  as  basis  functions  are  shown  as 


examples. 


Figures  2.22  and  2.23  show  the  solutions  for  a  10.5  X.  wide  PEC  strip  obtained  by  MoM  for 
the  normal  and  30  degrees  above  the  grazing  incident  cases,  respectively.  Notice  that,  since  MoM 
generates  a  fairly  good  solution  with  the  allowed  number  of  basis  functions,  the  improvement  in 
the  solution,  achieved  by  the  process  of  stepping,  is  not  as  noticeable  as  in  the  case  of  TM  strips 
(see  Figure  2.24).  Figures  2.25  and  2.26  show  the  results,  obtained  by  MoM  and  by  the  stepping 


AjnUCUiUP.M  ■ 


THETA  (DEGREES) 


Figure  2.22  Top:  Magnitude  of  the  induced  current  for  normal  incidence. 
Bottom:  Bistatic  RCS  for  normal  incidence. 

Solid  line  —  CGM. 

Dotted  line  --  MoM  using  20  plane-wave-type  functions. 
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Figure  2.23  Top:  Magnitude  of  the  induced  current  for  150  degree  incidence. 
Bottom:  Bistatic  RCS  for  150  degree  incidence. 

Solid  line  —  CGM. 

Dotted  line  —  MoM  using  20  plane-wave-type  functions. 


CURRENT 


W  Wil  W~W  WK  W\T 


X  POSITION 


31.52 


21.52 


co 

•a  11.52 
2 


-8.48 


-18.48 


THETA  (DEGREES) 

Figure  2.24  Top:  Magnitude  of  the  induced  current  for  30  degree  incidence. 
Bottom:  Bistatic  RCS  for  30  degree  incidence. 

Solid  line  ---  CGM. 

Dotted  line  —  Block  iteration  using  20  plane-wave-type  functions 
stepping  from  50  to  30  degree  incidence  with 
angular  increment  of  5  degrees. 
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procedure,  respectively,  for  a  10.5  X  wide  edge-loaded  strip  with  the  edge  load  width  1  X  from 
each  edge  of  the  strip  and  Rmax  equal  to  200  ohms.  The  improvement  achieved  by  stepping  is 
evident  from  the  figures. 


3.  CONCLUSIONS 


An  alternative  to  directly  inverting  the  large  MoM  matrix  is  to  recast  the  problem  into  a  form 
that  is  suitable  for  solution  via  iterative  schemes.  In  this  report,  a  number  of  iterative  techniques 
that  appear  to  be  useful  for  solving  the  problem  of  electromagnetic  scattering  from  electrically  large 
bodies  have  been  examined.  It  has  been  pointed  out  that,  while  the  use  of  iterative  methods  may 
enable  one  to  treat  larger  scatterers,  most  of  them  are  not  well-suited  for  handling  multiple 
excitations  in  an  efficient  manner.  Some  variational-iteration  schemes  based  on  the  use  of 
prechosen  entire  domain  basis  functions  that  are  suitable  not  only  for  treating  larger  bodies  but  for 
handling  multiple  incident  angles  as  well  have  been  suggested. 

When  compared  to  single-incidence  oriented  iterative  algorithms,  e.g.,  the  conjugate  gradient 
method,  variational-iteration  schemes  seem  to  converge  less  rapidly.  To  enhance  the  computational 
efficiency  of  the  variational-iteration  methods,  the  scheme  of  angle-stepping,  which  provides  better 
initial  guesses,  and  the  block-iteration  approach,  which  is  a  hybrid  method  based  on  a  combination 
of  iteration  with  Mom,  have  been  introduced.  Some  numerical  examples  have  been  given  to 
demonstrate  the  effectiveness  of  the  methods. 

Although  traveling-wave-type  functions  have  been  employed  as  basis  functions  in  this  report, 
the  methods  developed  here  can  be  applied  equally  well  with  any  type  of  basis  functions  that  can 
suitably  represent  the  induced  current  on  the  scatterer.  It  is  obvious  that,  for  efficient  computation, 
the  functional  form  of  the  prechosen  basis  functions  should  be  simple  so  that  they  can  be  easily 
implemented  in  the  computer. 
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